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The recently introduced two-parameter infinitely-many neutral alleles model 
extends the celebrated one-parameter version, related to Kingmans distribution, 
to difTusive two-parameter Poisson-Dirichlet frequencies. Here we investigate the 
dynamics driving the species heterogeneity underlying the two-parameter model. 
First we show that a suitable normalization of the number of species is driven 
by a critical continuous-state branching process with immigration. Secondly, we 
provide a finite-dimensional construction of the two-parameter model, obtained 
by means of a sequence of Feller diffusions of Wright-Fisher flavor which feature 
finitely-many types and inhomogeneous mutation rates. Both results provide 
insight into the mathematical properties and biological interpretation of the two- 
parameter model, showing that it is structurally different from the one-parameter 
case in that the frequencies dynamics are driven by state-dependent rather than 
constant quantities. 
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1 Introduction 

The two-parameter infinitely-many neutral alleles model is a family of infinite dimensional 
diffusion processes, introduced by Petrov (2009) and further investigated by Ruggiero and 
Walker (2009) and Feng and Sun (2010), that extends the celebrated one-parameter version, 
formulated by Watterson (1976) and characterized by Ethier and Kurtz (1981). Throughout 
the paper we will refer to the one- and two-parameter infinitely-many neutral alleles model 
simply as the one- and two-parameter model. More specifically, let 

(1) Voo = |z G [0, 1]°° : 21 > Z2 > . . . > 0, ^ < l| 

1=1 ^ 
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be the closure of the infinite-dimensional ordered simplex, and define, for constants < a < 1 
and > —a, the second-order differential operator 

i,j=l i=l 

where 6ij denotes Kronecker delta, acting on a certain dense subalgebra of the space C ( Vqo ) of 
continuous functions on Vqo • Then the closure of B generates a strongly continuous semigroup 
of contractions on C(Voo), and the sample paths of the associated process are almost surely 
continuous functions from [0,oo) to Vqo- Such process can be thought of as describing the 
temporal evolution of the decreasingly ordered allelic frequencies (zi, Z2, ■ ■ ■) at a particular 
locus in an ideally infinite population with infinitely- many possible types or species. See 
Feng (2010) for a review of infinite alleles models. Recently, Feng et al. (2011) determined 
the transition density for the two-parameter case. See also Borodin and Olshanski (2009) for 
a general construction related to Petrov (2009), and Ruggiero, Walker and Favaro (2013) for 
a partially related model with diffusive parameter 9. 

As shown by Petrov (2009) and Feng and Sun (2010), the two-parameter model is 
reversible and ergodic with respect to the Poisson-Dirichlet distribution with parameters 
{6, a), henceforth denoted PD(0,a). Introduced by Perman, Pitman and Yor (1992) (see 
also Pitman, 1995 and Pitman and Yor, 1997), this extends the one-parameter version 
PD{6) := PD(0, 0) due to Kingman (1975), and has found numerous applications in sev- 
eral fields. See for example Bertoin (2006) for fragmentation and coalescent theory. Pitman 
(2006) for excursion theory and combinatorics, Aoki (2008) for economics, Lijoi and Priinster 
(2009) for Bayesian inference, Teh and Jordan (2009) for machine learning and Feng (2010) 
for population genetics. 

Both these random discrete distributions arise as the decreasingly ordered weights of a 
Dirichlet process (Ferguson, 1973), when a = 0, and of a two-parameter Poisson-Dirichlet (or 
Pitman- Yor) process (Pitman, 1995) respectively. Alternatively, they can be constructed by 
means of the following so-called stick breaking procedure, also known as residual allocation 
model. Consider a sequence of random variables (Vi, V2, ■ ■ ■) obtained by setting 

n-l 

Vi = Wi, K = W„ J] (1 - W^), Wi Beta(l - a,d + m), 

i=l 

with < a < 1 and 9 > —a. The vector (Vi, V2, . . .) is said to have the GEM(^,q) distri- 
bution, named after Griffiths, Engen and McCloskey, while the vector of descending order 
statistics (V(i), V(2)i • • .) is said to have the PD(0, a) distribution. See Feng and Wang (2007) 
for an infinite-dimensional diffusion process related to GEM distributions. 
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Besides sharing the above stick-breaking construction strategy, it is weh known that the 
difference between these two random discrete distributions is structural and does not simply 
rely on a different parametrization. For example, the distribution PD(0) can be obtained 
by ranking and normalizing the jumps of a Gamma subordinator, whereas the PD(0, a) 
is obtained by performing the same operation on the jumps of a stable subordinator and 
appropriately mixing over the law of the normalizing factor (Pitman, 2003). See Section 2. 
Furthermore, the PD(0) distribution is obtained as weak limit of a Dirichlet distributed vector 
of frequencies (Kingman, 1975), while a similar construction for the two-parameter case is not 
available. For what concerns their diffusive counterparts, the properties of the one-parameter 
model, related to the PD(^) distribution, are well understood, whereas numerous are still the 
questions regarding the two-parameter model. In particular, given the above considerations, 
it is not surprising that a finite-dimensional construction of the process with operator (2), 
in terms of a sequence of finite dimensional diffusion processes, is currently available only 
when a = 0. To be more precise, consider the usual approximating diffusion for the Wright- 
Fisher discrete genetic model with n selectively neutral alleles and symmetric mutation. This 
corresponds to the operator 

i,j=l i=l 

acting on a suitable subspace of C^(\/n), with 

Vn = <Z e Voo : Zn+l =0, ^ Zj = 1 L 

i=i ^ 

with drift components 

(4) bP{z) = ^{l-z,)-ez„ 9>0. 

Ethier and Kurtz (1981) formalized the conditions under which the sequence of processes 
with operators defined by (3)-(4) converges in distribution to the one-parameter model, with 
operator obtained by setting a = in (2). As anticipated, a similar construction for the 
case < a < 1, 6 > —a is currently unavailable. It is to be said that two different sequen- 
tial constructions of the two-parameter model are given in Petrov (2009) and Ruggiero and 
Walker (2009). In Section 3.1 we will argue that despite offering interesting reads of the two 
parameter model, neither of these provides particular insight for the interpretation of the 
species dynamics underlying the infinite-dimensionality structure. In particular this is due 
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to the fact that both are based on finitely-many items. The problem at hand could then be 
rephrased as that of understanding from which Wright-Fisher-type mechanism, if any, the 
two-parameter model comes from. While the importance of providing a particle construction 
lies in the fact that the individual dynamics are dealt with explicitly, the contribution of 
a sequential construction by means of Wright-Fisher-type diffusions lies in the genetic in- 
terpretation one yields from the specification of the mutation rates at the nth step of the 
sequence. Such interpretation is clear in the case of (4), whereby each type has the same 
chance of mutating (cf. also (16) below), but is somewhat obscure for what regards the role of 
a in (2), especially in terms of its effect on finitely-many types. This role would be, at least 
partially, clarified by identifying suitable mutation rates, which are basic building blocks of 
the model and give important information on the reproductive mechanism of the underlying 
population. Historically, the (chrono) logical process has been the opposite, namely diffusion 
approximations were introduced for dealing mathematically in a simpler way with multi-type 
discrete models such as Wright-Fisher processes. But the recent advances, often stimulated 
by neighboring research fields, sometimes provide first what used to be the ultimate goal, 
as is the case with the two-parameter diffusion model, thus leaving a gap to be filled, both 
mathematically and in terms of interpretation. 

Motivated by these considerations, the purpose of this paper is to investigate what lies 
underneath the infinite-dimensionality of the two-parameter model in terms of the forces 
driving the species dynamics. We pursue this task in two different ways. First we derive 
an a-diversity diffusion for the two-parameter model. This is a continuous-time continuous- 
state extension of the corresponding notion for Poisson-Kingman models (Pitman, 2003), 
and describes the dynamics of the suitably normalized number of species in the underlying 
population. We show that such diffusion for the two-parameter model is a critical continuous- 
state branching process with immigration. This is done in Section 2, where a corresponding 
quantity for the one-parameter case is also discussed. The second goal is concerned with 
finding explicit transition rates for the mutation process which gives rise to the two-parameter 
model, and with providing a sequential construction for the limiting process in terms of finite- 
dimensional diffusions similar to those given by (3)- (4) for the one-parameter case. In Section 
3.1 we collect some brief considerations on the problem and on the fact that the existing 
constructions do not provide enough insight from a biological point of view. In Section 3.2 
we identify an instance of such mutation rates, which turn out to be unbounded on the 
usual domain of the process and to depend on the current species abundances. By means of 
some additional restrictions, we formalize a sequential construction where each term of the 
sequence is a Feller diffusion on a finite-dimensional subspace of Voo • 
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In achieving the two aforementioned goals, we are able to highlight a key difference be- 
tween the one- and two-parameter model, conveniently summarized by saying that the species 
dynamics of the former are driven by constant terms, whereas those of the latter are driven 
by density-dependent quantities. 



2 Heterogeneity in the two- parameter model 

The notion of a-diversity was introduced by Pitman (2003) for exchangeable partitions in- 
duced by random discrete distributions of Poisson-Kingman type. Let PK{g\t) denote the 
distribution of the weights {Pi)i>i = (Jj/T)j>i determined by the ranked points Jj of a Pois- 
son process with Levy density g, given T = t. A Poisson-Kingman distribution with Levy 
density g and mixing distribution 7 on (0,oo), denoted PK{g,^), is defined as the mixture 



/•oo 

pk(£',7) = / PK(^|t)7(dt). 
Jo 



For instance, the PD(0,a) distribution is obtained as a PK{ga,^e^a) model, for < a < 1 
and 6 > —a, where ga is the Levy density of a stable subordinator of index a, jg^a is 

and fait) is the density of a positive stable random variable of index a. Given an exchangeable 
random partition of N induced by a Poisson-Kingman distribution, i.e. such that its ranked 
class frequencies have distribution pk(^3, 7), this is said to have a-diversity S if and only if 
there exists a random variable S, with < 5* < 00 almost surely, such that 

(5) lim — = S a.s., 

where Kn is the number of classes of the partition restricted to {1, . . . ,n}. For instance, in 
the case of a PD(0, a) partition, we have S = T~" where T has distribution je^a- See Pitman 
(2003), Proposition 13. 

The idea of extending the concept of a-diversity from random distributions on simplices to 
a continuous-time continuous-state framework has been formulated in Ruggiero, Walker and 
Favaro (2013), where a certain rescaled, inhomogeneous random walk on the integers, which 
tracks the dynamics of the number of species in a normalized inverse-Gaussian population, is 
shown to converge to a certain one dimensional diffusion process on (0, 00). Here we derive 
an a-diversity diffusion for the two-parameter model, with the aim of providing insight into 
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the species dynamics underlying the infinite-dimensional process. The derivation is based on 
the particle construction given in Ruggiero and Walker (2009), here briefly recalled for ease 
of the reader. Let X^") = (Xi, . . . be a sample from a two-parameter Poisson-Dirichlet 
process, or equivalently (cf. Pitman, 1995) from the generalized Polya urn scheme given by 
Xi ~ Pq and 

(6) P{X,+i e-\X,,...,X,} = ^-^^Poi-) + ^ - a) 5x;{-) 

i=i 

for i = 2,...,n — 1. Here Pq is a non atomic probability measure on the space of the ob- 
servables (e.g. a Polish space), Ki < i denotes the number of distinct elements {X*, . . . , X'^J 
observed in {Xi, . . . ,Xi) and 5x* is a point mass at Xj. A simple way to make the sample 
into a Markov chain with fixed marginals is the following. Let X^"^ be updated at discrete 
times by replacing a uniformly chosen coordinate. Conditionally on Kn = k at the current 
state, and exploiting the exchangeability of the sample, the incoming particle will be of a 
new type with probability {6 + akr)/{9 + n — 1), and will be a copy of one still in the vector 
after the removal with probability (n — 1 — akr)/{9 + n — 1), where kr is the value of k after 
the removal. 

The following proposition recalls, in a discrete parameter version, the relation between 
the above described particle chain {X^"'\m),m € No} and the two-parameter model. For 
notational simplicity we omit here the details about the domain of the limiting operator 
(cf. (32)-(34) below). Here and throughout "=^" denotes convergence in distribution and 
Ca([0,oo)) denotes the space of continuous functions from [0,oo) to the space A. 

Proposition 2.1. [Ruggiero and Walker, 2009] Let Z(-) be the two-parameter model, cor- 
responding to the operator B as in (2). Let also {X(")(m),m G No} he the particle chain 
described above, and define y(")(-) = {y(")(t),t > 0} by y(")(t) = r]{X^''\[nH\)) , where 
rj{x^'^^) = (zi, . . . , Zn, 0, 0, . . .) if Zi is the relative size of the i-th largest cluster in x^"^ . Then 

y(-)(-)^Z(.), ^nC^^([0,oo)), 

as n —)■ oo. 

Hence the Markov chain {X^^\m),m € No}, once appropriately transformed and rescaled, 
provides a Moran-type particle construction of the two-parameter model. Denote now by 
{Kn{m),m G Nq} the chain which keeps track of the number of distinct types in X^'^\m), 
and let Mi^„ denote the number of types in X^'^\m) with only one representative. The 
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transition probabilities for Kn{m), denoted for short 



p{k, k') = F{Knim + 1) = k'\Kn{m) = k} 



are given by 



(7) 



p{k, k' 



1 



Min\ e + ak 



n 



+ n 



n{e + ny ^ " 

1 -p{k,k + 1) -p{k,k- 1), 
I 0, 



ifk' = k + l, 

iik' = k-l, 
if k' = k, 
else, 



for 1 < k < n. Here Mi^„/n is the probability of removing a cluster of size one, and k = 1 and 
k = n imply p(l,0) = and p{n,n + 1) = respectively. Since (7) need not be Markovian, 
we use an approximation of Mi^„ based on the following asymptotic result. From (5) and 
Lemma 3.11 in Pitman (2006), we have that the number Mi^n of clusters of size one observed 
in the sample is such that 

Ml, 



aS 



a.s.. 



so that Mi^n ~ Oik. This yields 

ak\ 9 + ak 
n J 6 + n 



(8) 



p{k,k') 



1 



+ o{n' 



ak 



-(n-l-a(A:-l)) + o(n-^+") 
1 - p{k,k + 1) -p{k,k- 1), 
0, 



iik' = k + l, 

ifk' = k-l, 

if k' = k, 
else. 



The following theorem identifies the a-diversity diffusion for the two-parameter model. De- 
note by Co([0,oo)) the space of continuous functions on [0,oo) vanishing at infinity. 

Theorem 2.2. Let {Kn{m),m £ Nq} be a Markov chain on N with transition probabilities 
as in (8), for < a < 1 and 6 > —a, and define {Kn{t),t > 0} by letting Kn{t) = 
Kn{\n^^°'t\) / n"' . Let also {Sg^a{t),t > 0} be a diffusion process on [0,oo) driven by the 
stochastic differential equation 



(9) 



dSe,a{i) = ^dt + \ 2aSe,a{t) dB{t), 
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where B{t) is a standard Brownian motion. If Kn{0) =^ ^^^^(O), then 
(10) {Kn{t),t>0}^{Se,a{t),t>0}, mqo,oo)([0,oo)), 



as n —)• oo. 



Proof. Denote by Un the semigroup induced by (8). For notational brevity, here we do not 
distinguish between (n, k) and (n — 1, /c — 1), since they are asymptotically equivalent. Then, 
for / G Co([0, oo)), we can write 



(Un - I)f 



f 



+ 



k + l 



/ 



/ 



k-l 



By means of a Taylor expansion we get 



(Un -I)f 



n 



a 9,a,k,n 



f 



k 



+ 



k 
1 



^ ak\9 + ak 
n J 6 + n 
ak{n — ah) 
niO + n) 



+ o(n-i+") 



where 



and 



(1) 

,a,k,n 



ak 
n 



2n2" 

+ ak ak / n — ak 



+ n 



n 



+ n 



+ o(n-^) 



(2) 



9,a,k,n 

Using (5), it follows that 
(11) 



^ ak\9 + ak akfn — ak 
n J 9 + n n \ 9 + n 



9 + n 
2ak 



9 + n 



+ o(n~^+"). 



sup |£/(s)-n^+°(CZ„-/)/(s)| ^0, as n ^ oo, 

sG[0,oo) 

where Cf{s) = 9f'{s) + asf"{s) is the infinitesimal operator corresponding to (9). Here (11) 
holds for every / belonging to an appropriate restriction Si{C) of Co([0, oo)) (to be formalized 
in Proposition 2.3 below). Under these conditions. Theorem 1.6.5 in Etliier and Kurtz (1986) 
implies that 



fl2) 



sup |C/(t)/(s)-t/„([ni+"tJ)/(s)|^0, /gCo([0,oo)), 

se[o,oo) 



as n — )• oo and for all t > 0, where U is the semigroup operator corresponding to C 
The assertion of the theorem with C[o ^^^([O, oo)) replaced by Z)[o_oo)([0) oo)) now follows 
from (12) and Theorem 4.2.6 in Ethier and Kurtz (1986). Finally, the convergence holds 
in C[o_oo)([0) oo)) since the limit probability measure puts mass one on C[o^oo)([0) oo)), and 
the Skorohod topology relativized to C[o^oo)([0) oo)) coincides with the uniform topology of 
C[o,oo)([0,oo)) (cf. Billingsley, 1968, Section 18). □ 
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Hence the dynamic heterogeneity of the two-parameter model is described by a non nega- 
tive diffusion obtained with a space-time rescahng which depends on the parameter a. Note 
that SQ^a{-) in (9) can be seen as a critical continuous-state branching process with immigra- 
tion (Kawazu and Watanabe, 1971; Li, 2006), obtained for example as diffusion approximation 
of a Galton- Watson branching process with immigration, with unitary mean number of off- 
spring per individual. Here 9 > is interpreted as the immigration rate, while the case 9 < 
has been treated in Going- Jaeschke and Yor (2003). 

The next proposition, which provides the complete boundary behaviour of the process 
driven by (9) and formalizes its well-definedness, is not new and included for formal com- 
pleteness. Let C be the second order differential operator 

d d2 

(13) C = 9— + as—,, 0<a<l,9>-a. 

as ds"' 

Proposition 2.3. The process {Sg^ait),t > 0} driven by (9) has the following boundary 
behavior: the boundary s = is absorbing for 9 < 0, instantaneously reflecting for < 9 < a, 
and entrance for 9 > a; the boundary s = oo is natural and non attracting for 9 < a, and 
natural and attracting for 9 > a. Moreover, Sg^ait) is null recurrent for 9 = a and transient 
for 9 a. For C as in (13), define 

&{C) = {/ e Co([0,oo)) nC2((0,oo)) : Cf G Co([0,oo))}, 



and 



/G^(£), if9>a, 

f G &{C) : lim,_,o x^/°/'(x) =0, ^/ < < a, 



/ G &{C) : lim^^o Cf{x) =0, if - a < 9 < 0. 

Then {{f,Cf) : f £ &Q^a{C)} generates a Feller semigroup on Co([0,oo)). 

Proof. The first assertion follows from Ikeda and Watanabe (1989), Example IV. 8. 2 and 
Karlin and Taylor (1981), Table 15.6.2. The second assertion follows from the first assertion 
and Ikeda and Watanabe (1989), Theorem VI. 3.1. The third assertion follows from the first 
assertion, together with the fact that exp{J-^'(0/a?/)dy} = and with Theorem 8.1.1 and 
Corollary 8.1.2 in Ethier and Kurtz (1986). □ 

The lack of positive recurrence immediately determines the non stationarity of the process 
{Se,a{t),t>Q}. 

We conclude the section with a brief discussion of the corresponding process for the one- 
parameter model. Although the notion of a-diversity is given for Poisson-Kingman models 
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with < a < 1 (cf. Pitman, 2003), a result analogous to Theorem 2.2 can be nonetheless 
derived for the one-parameter case, for which a = 0. The limit corresponding to (5) when 
a = is provided by Korwar and Hollander (1973) and is 



Kn 



lim 

n-!>oo log n 



a.s. 



Hence we expect the process for the normalized number of species to converge to a constant 
process, i.e. 



(14) 



{Kn{[cnt\)/\ogn,t> 0} Se,oit) 



for some c„ — )• cxd. Setting a = in (7) and proceeding similarly to the proof of Theorem 2.2, 
we get 



(Un - I)f 



k 



' -f( ' 



log n J log n \ log n 

w 



1 



w 



w n 



2 log^ n' \\ogn 
where w stands for the fact that 



1 



n J t) + n ^ n\0 + n 



n J + n n \ + n 
w / n 



-|-o((nlogn) ^] 



Mi^n ^ ~ Poisson(0); 

see Arratia, Barbour and Tavare (2003), Theorem 4.17. Hence in the limit the argument of 
the derivatives is constant, and Cn(Un — I)f{k/\ogn), with = nlogn, converges to 0. It 
follows that the dynamics of the number of species underlying the infinite-alleles models are 
driven by the constant process S'e^o(^) = ^ in the one-parameter case, and by the diffusion 
process So^ait) on [0, oo) with state-dependent volatility in the two-parameter case. This 
confirms the structural difference between one- and two-parameter models also from this 
dynamic viewpoint. A similar difference between the two cases will be found again in Section 
3.2 with a different approach. 



3 Finite-dimensional construction of the two- parameter model 
3.1 Preliminary remarks 

In the Introduction it was mentioned that two different sequential constructions of the two- 
parameter model have been provided in Petrov (2009) and Ruggiero and Walker (2009). In 



Species dynamics in the two-parameter diffusion model 



11 



this section we briefly outHne why these offer only partial insight into the dynamics underly- 
ing the two-parameter model from a biological perspective, motivating the need for further 
investigation. 

The above-mentioned constructions are given respectively by a sequence taking values in 
the space of partitions of N, and by the Moran-type particle representation outlined in Section 
2 above. Both cases are based on a dynamic system of finitely-many exchangeable particles 
and exhibit right-continuous sample paths, whereas (3) (with an appropriate domain) charac- 
terizes an n-dimensional diffusion process. Another notable feature of these constructions is 
the assumption that the mutant offspring distribution, namely the distribution that generates 
the mutant types, is non-atomic and thus selects types which appear for the first time with 
probability one (in the framework of Petrov (2009) this amounts to say that a new box is 
occupied with probability one). In particular, such feature turns out to be the key for proving 
the weak convergence of the sequences to the two-parameter model (see e.g. Ruggiero and 
Walker (2009) after Remark 3.1). Such assumption of non atomicity cannot be applied in a 
construction similar to (3)- (4), because the mass of the distribution must concentrate on the 
enumerated types, in order to keep the maximum amount of species constant in time. To 
be more precise about this aspect, note first that the drift coefficients in the Wright-Fisher 
operator (3) are determined as 

Here qij is the intensity of a mutation from type i to type j, and diagonal elements are qu = 
— Yl'^jki^iji SO that {qij)ij=i^...,n is a square matrix with nonnegative off-diagonal elements 
and row sums equal to zero. In general the mutation rate q^J^ (z) can be thought of as state- 
dependent, but in many interesting cases only the dependence on n is needed. The drift (4), 
for example, is obtained by taking parent-independent symmetric mutations with rates 

(16) '^ = ^^ 

whereby when a type i mutates, the new type will be any of the other n — 1 types with 
equal chances, and controls how often on average mutations occur. In this case then the 
mutant type is chosen with uniform probability, and the mutant offspring distribution is 
discretely supported. Such derivation of the one-parameter model can be extended to have 
non-symmetric mutation (see for example Ethier and Kurtz (1981), Theorem 3.4), but the 
difference is not relevant for our purposes. 

Hence the two existing constructions for the two-parameter model feature finitely-many 
objects, potentially of infinitely-many types, and a diffuse mutation offspring distribution. 
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while the desired construction should feature infinitely- many objects of finitely-many types 
and a discretely supported mutation offspring distribution. 

Prom a mathematical point of view, ideally we seek mutation rates ql.'^\z) yielding, 
through (15), the ith component limit drift term 

(17) 

and the b^"'\zys satisfy the boundary conditions 

(18) b'['\z)>0, ifzi = 0, b'['\z)<0, ifzi = l, 
for z G An and 

(19) A„ = |ze[0,l]": Zi>0, = l|. 

However, obtaining (17) and (18) jointly is clearly not possible, since in = the drift 
should be non negative for all n but strictly negative in the limit. Since condition (18) is 
crucial for the well-definedness of the nth term of the sequence, the alternative strategy will 
then be to relax (17) to the weaker condition 

(n)/ .df{z) unif , ^^f(z) 



(20) Y.t\^f-^^-^-E(^^^+-) g^. , 

i=l * i=l * 

as n — )• oo, for a sufficiently large set of functions f{z). Obtaining rates qij\z) which yields 
drift terms satisfying (20), together with some additional restrictions concerning the volatility 
and the state space of the process, will then suffice to provide the desired convergence. 

3.2 Sequential construction 

Let n > 2 throughout the section, and let A„ be as in (19). Consider a sequence of real 
numbers {e„}n6N satisfying 

(21) < e„ < ^ Vn, e„ = o(n-^), 

n 

and define the compact subspace of A„ given by 

An,e„ = |z e [0,1]" : z, > en,Y,z, = l| , 
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where z E ^n,s„ implies Zi G 1 — (n — 7^ for all i. See Figure 2 below. Consider 
the second order differential operator 

i,j=l i=l 

with domain 

(23) ^(^„) = {/: /gC2(A„,,J}, 
where 

c2(A„,,J = {/ g C(A„,,J : 3f e C2(M-), /|a„„„ = /} . 
The covariance components in (22) are specified to be 

(24) alf{z) = {zi - en){6ij{l - nSn) - {zj - e„)) 

_ J {zi- en){l - (n - l)e„ - Zi), if i = j, 

[ -{Zi - en){zj - En), ^^iy^j- 

These can be seen as Wright-Fisher-type covariance terms restricted to A^^^^, since they 
vanish at Zi = En and Zi = 1 — {n — l)En- Additionally, consider the state-dependent mutation 
rates 

(25) ^S"(-) = ;^ + i;;^Ii-»p{-('-.--")}l. 

Before providing some considerations on the form of q^j\z), note that (15) yields the drift 
components 

(26) ^(^ + 1)^ 

-a [1 - exp{-(zj - En)}] . 

Here the first two terms of b["'\z) equal 

a 

(1 - (n - l)e„ - Zi) - 9{zi - En), 



n — 1 



and Zi = 1 — {n — l)e„ implies Zj = En for all j / i. Using the last two observations in (26) 
shows that 

(27) b^\z)>0 iiz, = En, b^\z)<0 if = 1 - (n - l)e„, 
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SO that bi{z) satisfies (18) restricted to A„^£,^. 

In order to provide some intuition on (25), we have to consider separately the constant 
and the frequency-dependent term. The former attributes equal chances of mutation to all 
species regardless of their abundance, as in the one-parameter model. To evaluate the effect 
deviation from (16), recall now that the limit operator (2) acts on functions 
defined on (1), where the frequencies have been ordered. The same ordering operation will be 
done before taking the limit of An, with the formal appearance of the operator unchanged, 
so that it is correct to think in terms of ranked frequencies. In light of this, the term j at 
the numerator of the second term in qj:'j\z) can be interpreted as an approximate indication 
of the size of the frequency Zj, a greater j implying a lower Zj. Hence mutations from i to 
j occur more frequently if Zj is relatively low, implying a redistributive effect. This has to 
be interpreted as a conditional mechanism, related to the probability of directing the mass 
of the ith-type individual to some species j i, conditional on the fact that such individual 
mutates. In order to evaluate the unconditional chances of mutation of ith-type individuals, 
consider now the rescaled state dependent term in qlj\z), namely 

(28) z-i[l-exp{-(z, -£„)}]. 

Recall that the range of values of Zi is determined by n through e„ and grows to [0, 1] as n 
diverges, and note that the clear non monotonicity of the quantity in (28) is displayed on 
such range only for n large enough. Figure 1 plots (28) as function of Zi for n = 10, 10^, 10^ 
and En = n^^'^, so that (21) holds. The plot highlights the contribution of the rescaled state- 
dependent term of ql:j\z) with respect to the one-parameter mutation rate (16). The behavior 
of q^^^ (z) for Zi relatively far from e„ can be interpreted in terms of a reinforcement mechanism 
similar to that featured by the PD(0,a) distribution (see Lijoi, Mena and Priinster, 2007). 
It can indeed be observed that the probability that a further sample from (6) is an already 
observed species is not allocated proportionally to the current frequencies. The ratio of 
probabilities assigned by (6) to any pair of species is r^j- = (rij — a)/{nj — a). When 
a — )• 0, the probability of sampling species i is proportional to the absolute frequency n,, 
or equivalently to Zj, which in continuous time is reflected by a constant mutation rate as 
in (16). However, since rij is increasing in a, a value of a > reallocates some probability 
mass from type j to type i, so that, for example, for rij = 2 and nj = 1 we have rij = 2, 3, 5 
for a = 0, 0.5, 0.75 respectively. Thus a has a reinforcement effect on those species that have 
higher frequency. On the other hand, the behavior of q^J^ (z) for Zi near the boundary is what 
ultimately makes the process well-defined in a bounded region. For Zi | £n, (25) converges 



Species dynamics in the two-parameter diffusion model 



15 




Figure 1: Contribution of the (rescaled) state-dependent part of the mutation rate ^^"^(z) in 
(25) with respect to the rate of the one-parameter model. The plot shows (28) as 
function of Zi, with e„ = n~^'^ and, from innermost to outermost: ni — 10, n2 = 
10^ ns = 10^ 



to the rate of the one-parameter model (16), so that the associated drift behaves locahy as 
for the one-parameter case and Zi is kept inside the boundary e„ > 0. Roughly speaking, the 
decreasing part of the rate function pushes the frequencies towards smaller values, whereas 
the leftmost plotted part is responsible for keeping the frequencies inside the state space. 

The following result shows that the above defined operator characterizes a Feller diffusion 
on An^£„. Let || • || denote the supremum norm. 

Theorem 3.1. Let An be the operator defined by (22)-(23)-(24)-(26). Then the closure of An 
in C(A„^£^) is single-valued and generates a Feller semigroup {7^(t)} on C{/S.n,£n)- For each 
Vn £ =^(A„^e„), there exists a strong Markov process = > 0} corresponding 

to {7n{t)} with initial distribution z/„ and such that 

P{zH(.)gCa,_([0,oo))} = 1. 

Proof. It is easily seen that An satisfies the positive maximum principle on A^^^^, that is 
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if / G ^{An) and zq £ An,e„ are such that f{zo) = \\f\\ > 0, then Anf{zo) < 0. This 
is immediate in the interior of A„^e^, while on the boundaries it follows from (27) and the 
fact that (24) vanishes at every boundary point. Denote now z'^ = z^^ ■ ■ ■ z^" and a — 5i = 
{ai,...,ai-l,...,an) for ai, ct„ G N. Then 



Art.Z 



2^ n-1 n-1 n(n + l)^^ 



az" + ae^"e ^^z° 



+ {cJi - 1) {\\ -{n- \)enV'^^ - z'^ - e„[l - (n - \)enV~''^^ + Snz"-^^ 

Hence, letting denote the algebra of polynomials in (zi, . . . , z„) restricted to A„^e,j with 
degree not greater than m G N, the image of An computed on Lm contains functions belonging 
to Lm and of type e~^»z'^. Since for every g[x) G C{K), with K compact, and f{x) = e^g{x) G 
C{K), there exists a sequence {p^''^} of polynomials on K such that || / — p^^^^H — t- 0, so that 
||g-^p(*:) _ (^11 it follows that the image of An is dense in C(A„^e,J, and so is that of 

\ — An for all but at most countably many A > 0. Since UmLm is dense in C(An,£„), the Hille- 
Yosida Theorem (see Ethier and Kurtz, 1986, Theorem 4.2.2) now implies that the closure 
of An on C{An^e„) is single- Valued and generates a strongly continuous, positive, contraction 
semigroup {Tn{t)} on C(A„^e„). The fact that (1,0) belongs to the domain of An implies 
also that {Tn{t)} is conservative. Note now that for every zq G A„^£^ and 5 > there exists 
/ G ^{An) such that 

sup f{z) < f{zo) = 11/11 and A„/(zo) = 0, 

where B'^{zq, 6) is a ball of radius 6 centered at zq. Take for example f{z) = —Cs Y17=ii'^i ~ 
zq)^ for an appropriate constant Cs which depends on 6. Then the second assertion follows 
from Theorem 4.2.7 and Remark 4.2.10 in Ethier and Kurtz (1986). □ 

Given Vco as in (1), define the subspaces 

oo 

(29) Voo = G Voo : J2^' = l}' 

1=1 

and 

Vn,£„ = |z G Voo : Zn>en> Zn+l = o|. 
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Figure 2: Projection onto the first three coordinates of A„ (light blue), A„,e„ (dark blue), 
Vn,o (yellow), Vn,e„ (red), for n > 3. The (ordered) nth term of the sequence lives 
in the red region, so that the frequencies are bounded away from zero for all n. As 
n increases, the red region converges to the yellow region, so that the limit process 
is left free to move in the full ordered simplex. 



See Figure 2. Define also the Borel measurable map pn ■ A„ — )• Voo by 

(30) /9„(2;) = (z(i), . . . ,2;(„),0,0, . . .), z e An, 

where are the decreasing order statistics of z G A^. It is clear that p„ maps An^sn i'^to 
V„,e„. If Z(")(-) is the Markov process of Theorem 3.1, our aim is thus to show that 

(31) P„(ZW(.))^Z(.) 

in the sense of convergence in distribution in Cy ([0,oo)) as n — )• oo, where Z{-) is the 
diffusion process corresponding to the operator (2) with an appropriate domain. To this end, 
consider the symmetric polynomials 

(32) y^m{z) = J2 ^ e Voo, m>2, 

i>l 
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and define 

(33) ^{B) = |subalgebra of C(Voo) generated by 1, fsiz), if4{z), . . . |. 
Lemma 3.2. is dense in C(Voo)- 

Proof. In Ethier and Kurtz (1981) (see proof of Theorem 2.5) it is proved that the closure of 
^o('S) defined as 

(34) ^o{B) = jsubalgebra of C(Voo) generated by 1, (^2(-z), V23(-2), • • • | 
equals C(Voo)- Note now that 

zi= lim (/^^(z)^/™, Z2= lim (v9^(z) - zi)^/"", ... 

from which 



^2{z)= lim ^^{zfl-^ + {ip^{z)-z^fl^ + .^ 



SO that £ &{B). It follows that &q{B) = &{B), from which the result fohows. □ 

Before providing the convergence argument, we recall the relevant theorems about the 
formal existence and the sample path properties of the process Z(-) appearing in (31). 

Theorem 3.3. [Petrov, 2009] Let B be the operator (2) with domain (34)- The closure of 
B in C(Voo) generates a Feller semigroup {T{t)} on C(Voo); o-n-d for each v S =^(Voo) 
there exists a strong Markov process Z{-) = {Z(t), t >0} corresponding to {T{t)} with initial 
distribution v and such that 

P{Z(.)gC^^([0,oo))} = 1. 

Let Voo be as in (29). The following result shows that if the initial distribution of the 
Markov process Z{-) of Theorem 3.3 is PD(^, a), then the law of the process is concentrated 
onCv^([0,oo)). 

Theorem 3.4. [Feng and Sun, 2010] Let Z{-) = {Z(t),t>0} be the Markov process of 
Theorem 3.3, and assume Z{{)) ~ PD(0,a). Then 

P{Z(t) G Voo, Vt > 0} = 1. 



Species dynamics in the two-parameter diffusion model 



19 



Denote now by Bn the right hand side of (22), with coefficients as in (24) and (26) but 
domain 

(35) ^{Bn) = {fGC{Vn,eJ- / o p„ € ( A„,, J } . 

We are now ready to state the main result of the section. 

Theorem 3.5. Let Z^'^\-) and Z{-) be the Markov processes of Theorem 3.1 and Theorem 
3.3, with initial distribution Un G ^(A„^£,J and v G 9^{S oo)^ respectively. If ° Pn^ =^ i^, 
then (31) holds in Cy ([0, oo)). If in addition v is PD(0, a), then (31) holds in Cv^o ([0) co)). 

Proof. For p„ as in (30), define 7r„ : C(Voo) — ^ C{/S.n) by 7r„/ = f o and note that 
TTn : !^{B) ^^{An). Since for every / G ^{B) we have vr^iS/ = B{f o for all such 
functions and z G ^n,e„ we have 

AnTTnfiz) - TTnBfiz) = 



1 " 



dzidzj 



+ 



n 

EN 



z) + 6zi + a 



i=l 



dfiPnjz)) 

dzi 



It can be easily verified that the absolute value of the first term in square brackets is bounded 
above by a term of order 0(ne„), and that 



(36) \AnTTnf(z) - 7r„i3/(z)| < 0(ne„) J] 

df{Pn{z)) 



d^fiPniz)) 



dzidzj 



+ 0(n-^)E 



i=l 



dfiPniz)) 



dz 



i=l 



dz 



+ aEexp{-(zj - £„)} 



i=l 



dfiPniz)) 



dzi 



Observe now that for / G !^{B) of type (pmi x • • • x ipmki we have f{pn{z)) = f{z) and 

n k k n 



(37) 



E 

i=l 



df{z) 



dzi 



i=l j=l 



j=l i=l 



which is bounded above by Yl^j=i "^j < oo. Letting rrij = 3 (cf. (33)) the previous also implies 



Eexp{-(2;i - £„)} 



df(z) 
dzi 



uniformly as n — t- oo, where we have used (21) and the fact that the left hand side is maximized 
when Zi = £n. Furthermore, 
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is bounded, since the fact that Zi < i ^ makes the inner sum convergent as n — )• oo for 
ruj > 3. Finally 



(38) 



dzidzj 



rrihimh 



oo 



< 



+ '^mhmgiprnh-l<^mg-l 

q-^h i=/=h,q 

mh{mh - 1) + X] ^hrrig + ^ m,j 



whose right hand side is bounded. From (36), the above inequalities and (21) imply that 



(39) 



0, 



Given now that pniZ^"^^) clearly satisfies a compact containment condition (cf. Ethier and 
Kurtz, 1986, Remark 3.7.3), and that Theorem 3.3 implies that the closure ofB in Cy^ ([0, oo)) 
generates a strongly continuous contraction semigroup, the hypotheses of Corollary 4.8.7 in 
Ethier and Kurtz (1986) are satisfied and the first assertion of the Theorem follows with 
Cy^([0, oo)) replaced by Dy^([0,oo)), the space of cadldg functions from [0,oo) to Vqo- 
Furthermore, the convergence holds in ([0, oo)) C ([0, oo)), since the limit prob- 
ability measure is concentrated on Cy ([0,oo)) and the Skorohod topology relativized to 
([0, oo)) coincides with the topology on ([0, oo)). See for example Billingsley (1968), 
Section 18. Finally, the second assertion follows from Theorem 3.4 and from further rela- 
tivization to C^r^dO^oo)) of the topology on Cy ([0,oo)). □ 
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